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ABSTRACT 

The problem of radiowave propagation over irregular terrain is solved by using the 

wide angle parabolic equation. The terrain is characterized by its height profile and its 
ground constants (here conductivity a oo). We consider horizontal polarization and 
treat the ground as perfectly conducting (PEC) to simplify the formulation. This thesis 

uses a piece-wise conformal transformation to flatten the irregular terrain. The equations 
are solved by the split-step Fourier algorithm. A Hanning window is used in both the 
spatial and wavenumber domains to contain the computational domain. Effect of some 
numerical parameters such as the horizontal step size, and height of the computational 
domain on the accuracy of the solution is investigated. The numerical results are compared 
with available results for some typical propagation problems. 
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I. INTRODUCTION 



A. BACKGROUND 

Radiowave propagation over irregular terrain and in the presence of ducts is an 
extremely important topic in ground-to-ground as well as in ground-to-air 
communications. Similarly, the ability to predict radiowave propagation over irregular 
terrain has a significant impact in determining target detectability in a radar system. The 
physics of propagation, is affected by ever-changing atmospheric conditions and by 
complex terrain features on the ground. The path between the transmitter and receiver is 
often obstructed by natural or man-made obstacles, such as hills, buildings, atmospheric 
layers, rain, etc. In these cases waves can reach the receiver along more than one path and 
the phenomenon is termed as multipath. In the case of atmospheric multipath fading, 
super-refraction or sub-refraction can lead to abnormal radiowave propagation, which can 
result in high gain or severe loss of signal. Reflection multipath fading, which is due to 
interference between the direct and the ground reflected waves, depends strongly on the 
terrain geometry. It is important to everyone involved with communications and radar 
systems to predict the electromagnetic fields due to radiating sources in the troposphere 
taking into account ducting and terrain effects. 

Numerous analytical methods are available for predicting electromagnetic wave 
propagation, such as geometric optics, physical optics, normal mode analysis and 
combinations of the above. However, a complex environment complicates the application 
of some of these methods. The Parabolic Equation (PE) method has emerged as an 
extremely valuable method for assessing radiowave propagation in the lower atmosphere 
in the presence of ducts and over flat terrain, and for predicting radiowave propagation 
over sloping irregularities. One advantage of the parabolic Partial Differential Equation 
(PDE) is that the field at any location can be computed in terms of the field at a previous 
location. However, this would be accurate only when the waves propagate predominantly 
in the forward direction. This is approximately met in many propagation problems and is 
often assumed. 
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B. OBJECTIVE 

In this thesis we adopt a wide angle parabolic equation due to Thomson-Chapman 
[Ref. 1 ] to predict radiowave propagation over an irregular terrain and in the presence of 
ducts. We consider horizontal polarization and treat the ground as perfectly conducting 
(PEC) to simplify the formulation. In practical situations this assumption will have minimal 
effect for frequencies in the VHF band and above. The key feature of this thesis is to use 
piece-wise conformal transformations to map the irregular terrain into a flat one and use 
known techniques available for flat terrain. It is in the transformed space rather than in the 
physical space that the parabolic approximation made. Chapter II presents the derivation 
of the governing parabolic equation and the formulation of the conformal mapping. 
Chapter III details the generation of mesh in the transformed space and the numerical 
procedure for solving the PDE. The accuracy of the numerical solution is examined in 
Chapter IV. This includes a study on the effects of using different numerical values for 
various important parameters (e.g., step size, maximum height) on the accuracy of the 
solution and validation of the numerical results with known solutions for typical 
propagation problems. Recommendations and conclusions are presented in Chapter V. 
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II. FORMULATION 



In this chapter we present theory on the governing parabolic equation and 
conformal mapping. Figure 1 shows a horizontally polarized dipole placed over an 

irregular, PEC terrain. The terrain is characterized by its height profile and its ground 
constants (here the conductivity a ^ We wish to solve the fields at a point over the 
ground in the presence of irregularities and ducts. We consider a 2-dimensional case 

where the source, geometry, refractivity profile, and all fields are y-invariant. 



Tropospheric Botmdary z = z m a x 




Initial Range 



>x 




The terrain geometry is specified by means of discrete points, and straight line 
interpolation is assumed between the given data. The parabolic equation we present has to 
be finally solved numerically using either finite difference or Fourier transform approach. 
Both of these methods work best in a rectangular domain and a Cartesian mesh. Because 
the lower boundary is irregular, we cannot adopt a Cartesian mesh without discarding or 
adding points in going from one range step to the next. To avoid this situation, we will use 



3 



a coordinate transformation that will flatten the bottom boundary. Although there are 
multitudes of transformations available [Ref 2], we have decided to use conformal 
mapping in our case. The advantage of using a conformal mapping is that it preserves the 
Laplacian operator without introducing undesirable cross-terms. This will be particularly 
advantageous when split-step Fourier transform approach is used to solve the PE. In this 
thesis, we use the split-step Fourier transform approach to solve the parabolic equation 
due to its superior computational perfomance. Our development parallels that of [Ref 3]. 

A. WIDE ANGLE PARABOLIC PARTIAL DIFFERENTIAL EQUATION 

Figure 2 shows the coordinate system used in the thesis. 




Figure 2. Earth-centered spherical geometry. 

The starting point for our formulation is the earth flattened Helmholtz equation for 
the scaled field V [Ref 4]: 



Vl,V+kli\+2m)V=(i (1) 

where V(x,z)= sin 0 Ei^(x, z) for horizontal polarization, being the azimuthal 

component of the electric field, and (r,0,(t)) are the usual spherical coordinates. The center 

of the coordinate system is at the center of the earth and angle 0 is measured from the 
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vertical line joining the center and transmitter location. The transmitter is assumed to be 
located at (j)=0 (Figure 2). Furthermore, ko = (i>^ZoV^o , is the free space wavenumber, 
x=rsin0 is the range variable, z is the height variable, and m(x,z)=(n-l+z/aj is the modified 
refractive index of the atmosphere having an actual refractive index n, where a^ is the 
radius of the earth. In this chapter we include only the final details. Complete details can 
be found in [Ref 5]. 






A 



w=w(C) 





conformal mapping 



W- plane 






/-/■/ // / I 
^ - plane 



W=x + iz 



C = ^ + /Ti 



Figure 3. Transformation from W-plane to ^-plane. 



Under conformal mapping, [ Fig. 3], equation (1) is transformed to: 



V\^V+kl\w\Q\^(\+2m)V=0 



( 2 ) 



On expanding the Laplacian operator, equation (2) can be rewritten as 



— :r+ — r 



a^2 ari 



+kl\W {Q\^(\+2m)\=0 , 



( 3 ) 



where \W (Ql is the Jacobian of the transformation and is given by 
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(4) 



W(Q\ = l[ 



^+sin(§) 



V 



Substituting V=Ue'*°^, equation (3) can be simplified to: 



au 

-^=iko{Q-\)U 



(5) 



with Q-1- 




m , and we have made the reasonable assumption 



that I1F(Q|« 1. For the transformation considered in this thesis, this is true everywhere 
except near points P and Q of Figure 5. 



B. CONFORMAL MAPPING 



1. Conformal Mapping Approach 

A conformal mapping can transform an irregular polygonal shape into a 
rectangular one. We use a series of piece-wise conformal maps to transform the piecewise 
linear terrain to a flat one as shown in Figure 4. 




Figure 4. Terrain flattening. 
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For example, we conformally map the polygon AABB A shown on the left to a 
corresponding rectangle shown on the right. Because the slope is different from one 
segment to the other, the mapping is also different. Detailed derivations of the mapping 
are given in [Ref 5]. Here we present only the final results. 

2. Analytical Formulation of the Mapping 

Consider the transformation W=f(Q between the W=r+iz and ^=^+1^ complex 
planes as shown in Figure 5. 




Ar 
R' 



^=-Ar/2 

P' ^ Q' 

^ - plane 




IV=r + iz 



C = ^ + /ti 



Figure 5. Piece- wise conformal mapping. 



For the sake of simplicity we assume that the points R and S are at infinity. The 
mapping is given by the Schwarz-Christoffel formula [Ref 3]; 






^2 

!■ 



dk 



v+i i. 

0 X 2 ( 1->.)2 



( 6 ) 



where C 2 =^(sin(^)+ 1 ). 
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With reference to Figure 5, we define the slope prameter v as 



V = 



TC V/y/’ 




( 7 ) 



It is asumed that positive angles are measured in the counter clockwise direction. Points 
on the line PR are mapped to the points on the line PR through: 



zi = zj-i + ttJ(u;v) 



( 8 ) 



where as points on the line Q'S' are mapped to points on the line QS through: 



Ar 

Zr = J(ui~v) 



( 9 ) 



where 



t 2 y ^ A V , 

u=tanh (t^) =i — , and 

1-u’ 



u 

J(u;v)=l — 



dX 



0 X2^'’(l-X)2 



Accurate formulas for computing J(u;v) are given in [Ref 5], 



8 



III. SOLUTION PROCEDURE 



In the previous chapter we presented the underlying equations for the computation 
of fields over irregular/flat terrain and ducting conditions. In section A of the present 
chapter, we start with the generation of mesh points in the physical space. In section B, we 
continue with the numerical implementation using the split-step Fourier (SSF) algorithm 
for the final computation of the fields. 

A. GENERATION OF THE TRANSFORMED SPACE 

It is clear from equations (8) and (9) that a uniform mesh T|=kAr| in the 
transformed domain maps to different distribution of points on the vertical lines PR and 

QS. Furthermore, the distribution of points is dependent on the slope parameter v. 
Therefore, in going from one range step to another, the distribution of points on the right 
side of the previous step will not coincide with the distribution of points on the left side of 
the present step. This is illustrated in Figure 6 for v = 0.25, Ar| = 2 and Ar = 25 m, where 




Figure 6. Distribution of points about a slope discontinuity. 
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points just to the left of triangle base are indicated by circles and points to the right of it by 
inverted triangles. Clearly interpolation or extrapolation of data is needed before the field 
values calculated at one range can be used as input for the next step. In this thesis, we use 
rational approximation, [Ref 5] to do the desired interpolation. 

B. NUMERICAL IMPLEMENTATION USING THE SPLIT-STEP FOURIER 
(SSF) ALGORITHM 

Consider the parabolic partial differential equation given in equation (5). We 
would like to implement the equations by using the SSF algorithm. Figure 7 shows field U 
on a given vertical line. We would like to predict the field U+ on a subsequent line spaced 
Ar from the original line. 






u- 



777777777777777777 ^ V=0 

< > 



A r 



Figure 7. Marching in range. 

To solve the parabolic equation, the following Fourier transform pair is used in the 
transformed domain; 



Ufep) = f((7) = IufeT,)eW dT, 



(10a) 



(10b) 



—CO 



The boundary condition U(^,0)=0 for horizontal polarization implies that 

U(^,-ri)=-U(4,ri), ri>0, and U (^,-p) = -U (^,p), p>0. Using equation (10a) in equation 
(5) and making use of standard 'splitting of operators' [Ref 2] we have: 

U+ = U- ( 1 1 ) 



Arp^ 

where A= ^ . The field in the spatial domain is obtained by inverse 

ko+Jkl-p^ 

transformation and correcting for refractive effects: 

U+ = ' ( 12 ) 

In computing the ordinary Fourier transform pair, we use a N-point FFT. Let us 
assume that the various quantities are band limited over p^^ ^ p ^ Pm^x » that the 
transform is evaluated at p=0, Ap, 2Ap, ....(N-l)Ap where Ap=27t/(NAz). Positive 

wavenumbers occur at p=Ap, 2Ap, ,(N/2-l)Ap while negative wavenumbers occur at 

(N/2+l)Ap, (N/2+2)Ap, , (N-l)Ap. The value (N/2)Ap corresponds to both +p„„x 

■Pmax- *Jse a Hanning window to contain the computational domain at the top. The 
Hanning window is given by: 



h(n) = 



1 ,0<n<f 

^•2 ,4itns 'iN ^ ^ ^ N 

Sin (— ) >T-"-T 



(13) 



We use it both in the spatial and wavenumber domains. A mirror image of h(n) about n=0 
is used for negative wavenumbers. Note that the Hanning window forces a gradual rolloff 
to zero over the last quarter of the domain. With the use of the Hanning window the 
solution of the parabolic equation can be written: 

t/+ = [h(n)e-‘^Flh(n)U-]] (14) 
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As we mentioned at the begining of this chapter, because consecutive steps are 
taken in different transformed spaces, the calculated field at the right side of one range 
step cannot be directly used, in general, as input to the next step. For example, at the 
junction between flat and sloping terrain, points are uniformly distributed over a flat 
surface, but unevenly distributed over the sloping terrain. As mentioned previously, we use 
the rational approximation to calculate the field f(z) at a desired point z, in terms of the 
fields f ,f^. and f^ at positions z. , z^. and z^ respectively as shown in Figure 8. 



z 

z 



+ o 
o** 






z o 



Figure 8. Field interpolation. 



The field is interpolated as: 



= +ai(z-Zo))l(\^bi{z-Zo)). ( 15 ) 

where a i = (/"+5_ -f.s^)l(f+ -/_) , 6 1 = ( 5 = - s +)!(/+ -/_) and 5_ = (fo -/-)/(^o - z ) , 



•J+ = (f+-fo)l(z+-Zo). 



12 



IV. RESULTS 



Having introduced the basic theoretical background, we now present the numerical 
results using the developed algorithm for several practical cases of radiowave propagation. 
We start with the prediction of propagation factor (excess loss over the free space case) 
over flat earth, with and without a duct. We compare our results with the exact solutions 
or with results available in the literature. The next case we consider is that of propagation 
over spherical earth in a non-refractive atmosphere. Finally, we present our solution for 
propagation over a knife-edge shaped obstacle. We use this terrain irregularity to 
comment on the overall accuracy of the method and its dependence on parameters, such as 
the step size for the marching and the maximum height chosen to truncate the 
computational domain. 

A. PROPAGATION OVER A PEC PLANE 

Consider a transmitter and a receiver at a height of 5 meters above flat ground and 
separated by 1 km , as shown in Figure 9. 




Figure 9. Propagation between two stations A and B, both of which are stationed over 
perfectly conducting ground. 
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In our numerical approach, the propagation factor, PF, is obtained from the 
normalized field, V, as [Ref 6] 

PF=101og(lFl^xXo) (16) 

where x is the horizontal from the transmitter. A positive (negative) value of propagation 
factor implies gain (loss) with respect to propagation in free space. The exact solution can 
be obtained from image theory [Ref 7, p. 144], The expression valid in the far-field and 
for X » h,H is PF = \E! Eo\ = 2 s\n{k oh cosQ), where cos Q=H/-JlF +x'^ , h is the 
transmitter height, and H is the receiver height. The transmitter placed at x=0, z=5, where 

X, z are in meters. The frequency of operation is 1 GHz. We choose N=1024 , ziz=Ar|= 
=0.15 m and Ax=Ar=5 m. The height of the upper boundary is z^^,^ =76.5 m. The electric 
field is determined from the source line up to a distance of 1 km. 

The variation of the propagation factor versus the horizontal distance is plotted in 
Figure 10. It is seen that there is a very good agreement with the image theory solution. 
For the first 50 m, the numerical solution differs slightly from the known solution. This is 
due to the breakdown of the condition x » h, H in the image theory, far-field formula. 

B. PROPAGATION OVER A PEC PLANE IN THE PRESENCE OF DUCT 

Next we present the numerical results for radiowave propagation over a PEC 
plane, in the presence of a tri-linear duct [Figure 11] specified by 



’ 340-1-0.118Z 



M(z) = < 



499.03 -1.06z 
322.33 + 0.118Z 



0<z< 135 
135 <z< 150 
z > 150 



( 17 ) 



where z is in meters and M = m x 10®. This is a surface elevated duct with AM=16 and 
Ah=150 m, where AM and Ah are defined at Figure 11. The minimum frequency (cut off 
frequency) required for propagation in the presence of this duct is 217 MHz. In our case 
the frequency of operation of 3 GHz is well above this cut off frequency. The height of the 
transmitter is 30 m and the maximum height is chosen to be z =512 m. We choose 



14 



Propagation Factor (dB) 




Figure 10, Propagation factor versus the horizontal distance for a source, placed at a height of 
5 m. Flat earth, f =1 GHz , Ax=5 m and H=h=5 m. 
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height (m) 




modified refractive index M 



Figure 1 1 .Tri-linear surface elevated duct. 
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Receiver Height (m) 




Figure 1 2. Propagation factor versus the receiver height at a range of 40 km in the presence of a tri-linear 
duct. Flat earth, transmitter height =30 m, f =3 GHz. 
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N=512, Dz=2 m and Dx=200 m. Figure 12 shows propagation factor versus the receiver 
height at a distance of 40 km. A favorable agreement with [Ref 6, Figure 3] is obtained. 

C. PROPAGATION OVER SPHERICAL EARTH 

Figure 13 shows a transmitter and a receiver located beyond the horizon over 
spherical earth. We consider transmitter and receiver antennas at heights of 10 m each. 
Distance to the horizon, d.j. , from each station is JlhRe , where is the effective radius 
of the earth equal to 6375 km. In our case dj. =1 1.3 km. The situation described above is 

equivalent to the one where radiowaves propagate over flat earth but in the presence a 
modified atmosphere described by modified index 

m(z) = (n - I + z/ae) (18) 

where Oe is the radius of the earth, and we choose the refractive index, n, equal to 1 . We 
use the equivalent model in our numerical calculations. 



n=l 



m=z/ag 



transmitter receiver 





transmitter receiver 



1 


f Y 

h=10m H=10m 




^ 3^ 



d= 10 km 

Equivalent Problem 



Figure 13. Propagation over spherical earth. 

We consider a frequency of 300 MHz and treat the earth as a PEC. The exact solution to 
the problem can be obtained by using the mode theory [Ref. 8]. Our case pertains to the 
case of low antennas where the heights are less than the critical height h^ , defined by 
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(19) 



he = 30A,o in m . 

In our case h^=30 m at f=300 MHz. The distance d, must satisfy; X.„d > 2 h H, meaning 
that the two stations must be at least 200 m apart. We use N=1204, Ax=50 m and Az=0.5 
m. This, places the maximum height at z^3,^=255.5 m. We solve for the propagation factor 
as a function of the horizontal distance using our model and compare the result with the 
exact solution [Ref 8] 



PF =l£/£ol=272<l|^|^/„(/»)/-„(i?)l (20) 

2 

where 5 = [A:o/?e]3(src- 1), Sj^=s^-io, where is the relative permittivity and is the 
relative conductivity, Q = sd, s= {kolR\\ ^ = complex numbers which depend on 5, and 

4 (h), 4(H) are height gain functions for the n-th mode (term). For low antennas 4 (x) is 
independent of n and is given by /(x) = 1 +ikoxJsrc— 1 • The contants t,, are given as in 

Table 1: 



n 




1 


1.856e*”^ 


2 


3.245e'"'^ 


3 


4.382e‘”'^ 


>3 


0.5[37i(n+0.25)]^e‘'”^ 



Table 1 . Mode numbers for 5=oo ,(PEC). 

In Figure 14 we plot the propagation factor versus the horizontal distance d for a 
range up to 10 km. There is a very good agreement with the modal solution, beyond the 
first 400 m. It was already mentioned that the distance between the two stations has to be 
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Propagation Factor (dB) 




Figure 14. Propagation factor versus horizontal distance over spherical earth, f =300 MHz, Ax=50 m. 



greater than 200 m in the modal solution, which explains the disagreement for the first few 
hundred meters. 

D. PROPAGATION OVER PEC KNIFE-EDGE 

The final example we consider is that of propagation over a knife edge on a 
conducting plane. Consider a perfectly conducting knife-edge of height 25 m, as shown in 
Figure 15. 




Figure 15. Perfectly conducting knife edge between the transmitter at A and the receiver 
B, both of which are on perfectly conducting ground. 



We compare our numerical solution with a four ray model of knife edge diffraction 
theory [Ref 9]. The two stations are at distances d, and from the plane of the edge. We 
choose d,=d 2 =500 m. In the parabolic equation computations, we replace the knife edge 
with a triangular hill having a base of 50 m (Figure 15). The source is placed at a height 
of 3 m and we would like to determine the propagation factor as a function of receiver 
height. We choose a maximum height of 255.5 m, and the other parameters as: N=1024, 
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Az=0.5 and Ax=l m. The frequency of operation is 300 MHz. The numerical results are 
plotted in Figure 16, along with the four-ray solution, [Ref 9], There is an excellent 
agreement between the two plots. 

E. PARAMETRIC STUDY OF THE SSF ALGORITHM USING CONFORMAL 
MAPPING 

We use the knife edge irregularity to show the dependence of the solution on the 
following parameters; the horizontal step size and the maximum height chosen. We choose 
three different heights, each of which corresponds to the following values of N (with fixed 
vertical step size Az=0.5 m): N=256, 512 and 1024. The values of height are: z,„^=63.5, 
127.5 and 255.5 meters respectively. For each of these heights, we create three sets of 
plots, for three different values of the horizontal step sizes: Ax= 5, 2.5 and 0.75 meters.In 
Figures 17, 20 and 23 (Ax=5 m), we discuss the effect of the maximum height chosen for 
fixed Ax and it is seen that =76.5 m (Figure 17) is certainly not high enough. Only if 
we use higher a computational domain (Figures 20 and 23) do we obtain a solution which 
has resemblance with the four-ray solution. In the second set of Figures for a fixed Ax=2.5 
m, (Figures 18, 21 and 24), again we conclude that the higher the z^^ the better the 
agreement between the numerical and the four ray solution. In Figures 21 and 24 (z,^ 
equal to 125.5 and 255.5, respectively) we start to obtain a solution similar to the four-ray 
theory solution, but not yet an accurate one. Finally in Figures 19, 22 and 25, for the fixed 
horizontal step size Ax=0.75 m, we see that the solution with z^^ of 125.5 m or 255.5 m 
(Figures 22 and 25) produces very good agreement. 

Next we discuss the effect of the horizontal step size Ax for a fixed vertical step 
size Az and z^^ (N). In Figures 17 through 19 the maximum height is 76.5 m and Az=0.5 
m. We see that even for very small Ax, an acceptable solution is not obtained. In this case 
the solution is dominated by the height parameter (not convergent yet). In Figures 20 
through 22, we use z„, 3 ,,= 125.5 m and Az=0.5 m. It is seen that as we decrease the 
horizontal step size a better agreement with the four-ray solution is achieved (Figure 22). 
Figures 23 through 25 present the results for a fixed Az=0.5 m and z^^ =255.5 m. Again 



22 




Propagation Factor (dB) 



Figure 16. Propagation factor versus the receiver height at a range of 1 km. Smooth earth, the 
source placed at a distance of 500 m from the edge. Az=0.5 m, Ax=0.75 m and =255.5 m. 
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Receiver Height (m) 




Propagation Factor (dB) 



Figure 1 7. Propagation factor, PF, versus the receiver height. Az= 0.5 m, Ax=5 m and =63.5 m. 




Figure 18 Propagation factor versus the receiver height. Az= 0.5 m, Ax=2.5 m and =63.5 m. 
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Figure 19. Propagation factor versus the receiver height. Az= 0.5 m, Ax=0.75 m and =63.5 m. 
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Figure 20. Propagation factor versus the receiver height. Az= 0.5 m, Ax=5 m and =127.5 m. 
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Receiver Height (m) 




Propagation Factor (dB) 



Figure 21. Propagation factor versus the receiver height. Az= 0.5 m, Ax=2.5 m and =127.5 m. 



Receiver Height (m) 




Figure 22. Propagation factor versus the receiver height. Az= 0.5 m, Ax=0.75 m and -127.5 m. 
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Receiver Height (m) 




Figure 23. Propagation factor versus the receiver height. Az= 0.5 m, Ax=5 m and =255.5 m. 
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Figure 24. Propagation factor versus the receiver height. Az= 0.5 m, Ax=2.5 m and =255.5 m. 
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Figure 25 Propagation factor versus the receiver height. Az= 0.5 m, Ax=0.75 m and =255.5 m. 
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the effect of Ax is easily seen and Figure 25 shows a numerical solution with an excellent 
agreement with the four-ray theory solution. 
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V. CONCLUSIONS 



In this thesis, a numerically efficient method to model radiowave propagation over 
irregular terrain and in a stratified atmosphere, using the parabolic equation, was 
implemented and tested. The parabolic equation method is a full-wave method, and 
aspects such as forward reflection, refraction, diffraction, and surface wave propagation 
are included. However, it ignores back-scattering. Since it ignores back propagation, it 
allows for a rapid solution of the fields using marching techniques along the propagation 
path starting from the source. It offers an advantage compared to the ray method in that it 
is valid in the shadow region where the latter method completely breaks down. 
Furthermore, it is the most practical method for predicting propagation over long ranges 
(thousands of wavelengths). 

After presenting the wide angle parabolic equation, we presented a piece-wise 
conformal transformation to map the irregular terrain into a flat one, so that known 
techniques available for flat terrain can be used. The split-step Fourier (SSF) algorithm 
was used in the computational domain to solve the parabolic PDE. A Hanning window 
was used both in the spatial and wavenumber domains. Although very small step sizes 
were used for marching along ranges of several kilometers, the method is very time 
efficient. Horizontal polarization was considered and the ground was treated as a perfect 
electric conductor (PEC) to simplify the formulations. This assumption will have minimal 
effect in practical situations and for frequencies in the VHF band and above. 

Numerical results to predict radiowave propagation for several practical cases 
were computed and validated using the results available in the literature. These included 
both propagation over flat terrain in a refractive and a non-refractive atmosphere, 
propagation over spherical earth, and propagation over a PEC knife edge placed on a 
perfectly conducting plane. Excellent agreement was observed and demonstrated for all 
cases. Also included, was a parametric study of the algorithm used, to show the effect of 
the horizontal step size Ax and the maximum height ,for a fixed vertical step size Az. 
The need for larger heights and smaller horizontal steps was demonstrated. 

The model presented in this thesis can be applied to communication problems or 
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for target detection where prediction of propagation of HFA/HF signals over 
inhomogeneous and irregular terrain is highly desirable. Horizontal polarization applies 
mainly for radar systems. For this case, signals received from a target depend on the direct 
wave and the ground reflected wave. The latter ray's path depends on the terrain 
roughness and ground constants. It is also important to assess the atmospheric conditions. 
The numerical model developed in this thesis provides a useful tool for predicting the 
performance of a radar or communication system before its final development. 
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